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Several approaches to photonuclear reactions, based on the time-dependent 
density-functional theory, have been developed recently. The standard lin- 
earization leads to the random-phase approximation (RPA) or the quasiparticle- 
random-phase approximation (QRPA). We have developed a parallelized 
QRPA computer program for axially deformed nuclei. We also present 
a feasible approach to the (Q)RPA calculation, that is the finite ampli- 
tude method (FAM). We show results of photoabsorption cross sections 
for deformed nuclei using the QRPA and FAM calculations. Finally, the 
canonical-basis approach to the time-dependent Hartree-Fock-Bogoliubov 
method is presented, to demonstrate its feasibility and usefulness. 

PACS numbers: 21.10.Re; 21.60.Jz; 24.30. Cz 

1. Introduction 

Photonuclear reaction cross sections are the fundamental properties in 
nuclear systems. In the energy region of giant resonances (E = 10 ~ 30 
MeV), the absorption process is dominated by the electric dipole excita- 
tions. The giant dipole resonance (GDR) has been of significant interest 
in studies of nuclear structure and reaction. It exhausts almost 100 % of 
the energy-weighted sum-rule value, corresponding to a collective oscilla- 
tion of neutrons against the protons. A typical measurement of the GDR 
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in stable nuclei is the photoneutron cross section measurement using mo- 
noenergetic photons [!]. The energy of the GDR peak was found to have 
a mass dependence midway between A" 1 ^ and A~ l l & which correspond to 
the Steinwedel- Jensen and Goldhaber- Teller models, respectively [2]. The 
energy-weighted sum-rule value is larger than the classical Thomas- Reiche- 
Kuhn (TRK) value by 20 % pQ in average. The general trend of the width 
of the GDR is well correlated with the neutron magic numbers, which may 
suggest that the main origin of the spreading width is the shape fluctuations 
in the ground state H] • The double-peak structure in GDR appears for 
axially deformed nuclei, known as the deformation splitting, because of the 
different frequencies for vibrations along and perpendicular to the symmetry 
axis pp. 

In a microscopic point of view, one can construct the giant resonance 
from a superposition of particle-hole excitations. Since dynamics of the giant 
resonances are basically in a small-amplitude regime, the random-phase 
approximation (RPA) [2] has been extensively utilized for studies of their 
properties. Although the spreading width is not taken into account in 
the RPA level, main features of the giant resonance are well reproduced. In 
this paper, we present three theoretical approaches to studies of the nuclear 
response; the standard quasiparticle RPA (QRPA) [2j, the finite amplitude 
method (FAM) [5], and the canonical-basis time-dependent Hartree-Fock- 
Bogoliubov (Cb-TDHFB) method [UJ. The numerical results will be shown, 
mainly focused on the photoabsorption cross section. 



2. Quasiparticle random-phase approximation for axially 

deformed nuclei 

The quasiparticle RPA (QRPA) is a standard method to calculate lin- 
ear response in heavy open-shell nuclei [2]. However, since its application 
to realistic energy functionals requires a large computational task and a 
complicated programing, the QRPA calculation for heavy deformed nuclei 
is still a challenging subject at present. 

We have recently developed a parallelized computer code of the QRPA 
based on the Hartree-Fock-Bogoliubov (HFB) state with the Skyrme func- 
tionals, which is an extended version of that developed in Ref. [7J, to in- 
clude the residual spin-orbit interaction. A missing part is only the residual 
Coulomb interaction that does not significantly affect nuclear response func- 
tions (See Sec. I4.2p . 

First, we solve the following self-consistent HFB equation for the quasi- 
particle states: 

-A* A -(/* - A)*) (k) = E * UU (1) 
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Fig. 1. Calculated (lines) and experimental (symbols) pliotoabsorption cross section 
in 154 Sm. The SkM* parameter set and the smoothing parameter of T = 2 MeV is 
used. See text for details. Experimental values are taken from Ref. [3J. 



where the single-particle Hamiltonian h[p, k] and the pair potential A[p, k] 
are functionals of the density p and the pairing tensor k. The self-consistent 
solution of Eq. (P) determines the ground-state densities (pq,kq) and the 
ground-state Hamiltonians (/io,Aq). To describe the nuclear deformation 
and the pairing correlations, simultaneously, in good account of the contin- 
uum, we solve the HFB equations in the cylindrical coordinate space. We 
assume axial and reflection symmetries in the ground state. To reduce the 
QRPA-matrix dimension, we introduce a cut-off energy £c qp = 60 MeV for 
the two-quasiparticle states. For instance, the number of two-quasiparticle 
states becomes about 38,000 for the K w = _ excitation in 154 Sm. Then, we 
calculate the QRPA matrix elements and diagonalize the matrix, to obtain 
the QRPA normal modes. 

y ( A Pn& B Pn& \(X lS \ =f ^( X lS \ 

Since the spreading effect is missing in this calculation, the dipole strength of 
each discrete eigenmode is folded by the Lorentzian curve with a smoothing 
parameter T. A more detailed description can be found in Ref. |S]. 

We show in Fig.Q]the photoabsorption cross section for 154 Sm. The HFB 
calculation with the SkM* parameters produces the ground state in a prolate 
deformation of j3 = 0.31. It clearly shows a deformation splitting due to a 
prolate deformation of the ground state. The experimental data [S] are well 
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reproduced in the calculation. We have carried out a systematic analysis on 
Nd and Sm isotopes and have found that the spreading effect with T = 2 
MeV can well reproduce experimental data from spherical, transitional, to 
deformed nuclei [8] . Especially, the agreement on the evolution of the GDR 
width as a function of the neutron number is excellent. 

3. Finite amplitude method 

In this section, we recapitulate the methodology of the finite amplitude 
method (FAM) we have developed for small-amplitude oscillations based on 
the time-dependent density- functional theory [5j |9] . 

3.1. FAM without pairing correlations 

First, we discuss the case that the energy density functional is repre- 
sented by normal density p only. In this case, the density can be expressed 
by the Kohn-Sham orbitals, p = ^ \4>i)((fii\, where the subscript i indicates 
the occupied orbitals (i = 1,2,- •• , ^4). The linear-response equation to 
a weak external field with a fixed frequency, V ex t(u), can be expressed in 
terms of the forward and backward amplitudes, \Xi(uj)) and (Yi(u))\. 

to \Xi(u)) = (h - ei) \Xi(u>)) + P {V cxt (ui) + 5h(u>)} \<j)i), (3) 
-uj{Yi{u)\ = (Yi(a;)|(/io-ei) + (^|RxtM+^M}-P- (4) 

where the operator P denotes the projector onto the particles space, P = 
1 — Yli Usually, the residual field 6h(u) is expanded to the first 

order with respect to \Xi(u))) and \Yi(u)). This leads to the well-known ma- 
trix form of the linear-response equation, known as the RPA. For deformed 
nuclei, the calculation of these matrix elements is time-consuming in prac- 
tice and their storage requires a large memory capacity. In the FAM, we do 
not explicitly linearize the equations. Instead, we utilize the fact that the 
linearization can be numerically achieved for 5h(uo) = h[po + Sp(u))] — ho, 
if the transition density 5p(uj) is small enough to validate the linear ap- 
proximation. The FAM is nothing but a trick to perform this numerical 
differentiation in the single-particle (Kohn-Sham) Hamiltonian h[p\. 

The residual field 5h(uj) depends only on the forward "ket" amplitudes 
\Xi(oj)) and backward "bra" ones (Yi(u)\. In other words, it is indepen- 
dent of bras (Xi(ui)\ and kets \Yi(u))). This is related to the fact that the 
transition density Sp(u) depends only on \Xiioj)) and (Yi(u))\. 



5p(u) = ^ {|XiH>(&| + \<f>i)(Yi(u)\} 

i 



(5) 



draft printed on January 18, 2011 



■5 



Ca 










1 












XX4IXXXXXX 




Ar 




■ 


























.1 


X 


















X s 




' 






i 














XXI Si 


























Mg 




















Ne 






i 








i 













... 




, ... 


c 




— 




iv 




Be 








He 





Fig. 2. Electric dipole strength distribution in light nuclei. The horizontal axis 
corresponds to excitation energy of ~ 35 MeV. The SkM* parameter set and the 
smoothing parameter of T = 1 MeV is used. 



We calculate the residual field by introducing a small real parameter r/ to 
realize the linear approximation [5]. 

8h(u) = -(h [p n ] - h ) , (6) 

where ho is the Hamiltonian for the ground state and are defined by 

pv = E {(i&> + v\Mu)))m\ + ^hi)} . (7) 

i 

Once \Xi(u))) and (li(a;)| are given, the calculation of h[p v ] is an easy task. 
This does not require complicated programming, but only needs a small 
modification in the calculation of h[p\. Of course, eventually, we need to 
solve Eqs. ([3]) and dU to determine the forward and backward amplitudes. 
We use an iterative algorithm to solve this problem. Namely, we start from 
initial amplitudes |X! ) and (K |, then update them in every iteration, 

(\x\ n } ), (y. (n) |) (|X 4 (n+1) ), (i; (n+1) |), until the convergence. In each step, 
we calculate 8h{u) using the FAM as Eq. ([6]). For more details, readers are 
referred to the reference [5]. 

We have developed a parallelized computer program of the FAM for 
a Skyrme functional in the three-dimensional (3D) coordinate-space rep- 
resentation [9]. Currently, we are performing a systematic calculation of 
the electric dipole response in even-even nuclei. So far, we have calculated 
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Fig. 3. Calculated isoscalar monopole strength distribution for 50 Ca. The SkM : 
parameter set and the smoothing parameter of T — 1 MeV is used. 



the photoabsorption cross section in nuclei with A < 100. In Fig. [21 we 
demonstrate a part of our achievement for nuclei up to Ca isotopes. 

For nuclei with A < 40, the observed strength up to 30 MeV exhausts 
only 60 ~ 100 % of the TRK sum-rule value [TOj. This indicates that the 
considerable amount of the GDR strength is located above 30 MeV in light 
nuclei. We also observe that, although the RPA (FAM) calculation repro- 
duces a gross feature of the dipole strength distribution, it systematically 
underestimates the GDR peak energy by a few MeV for light nuclei [9]. 

3.2. FAM with pairing correlations 

The FAM in the previous section can be extended to superfluid nuclei, 
namely, to the QRPA with the HFB formalism. A self-consistent solution 
of Eq. ([1]) determines the ground-state densities (po,Ko) and the ground- 
state Hamiltonians (ho, Aq). Then, following the same argument as that in 
Ref. [5], we can derive equations for the residual fields, 5h(uj) and 5A(u) as 
follows: 

6h(uj) = - (h[p v , k„] - h ) , (8) 
V 

6A(u) = -(A[p v ,K v ]-A ), (9) 

where the density and pairing tensor (p rj , defined by 

p v = (V* +rjUX)(V + r ] U*Y) T , (10) 
Kr, = (V* +r,UX){U + r,V*Y) T . (11) 

Here, the forward and backward amplitudes (Xuv, Y^v) have subscripts pv to 
specify two-quasiparticles. On the other hand, the subscripts of (Uk^jV^) 
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indicate a basis of the single-particle space (k) and the quasiparticle (//). 
Again, utilizing an iterative algorithm for solution of the QRPA equation, we 
can solve the QRPA linear-response equation without explicitly calculating 
the residual interactions. 

We show in Fig. [3] an example of our FAM calculation for isoscalar 
monopole response in 50 Ca. We use the same parameter set and the same 
pairing energy functional as those in Ref. [llj. The quasiparticle states are 
truncated by the maximum quasiparticle energy of E qp = 200 MeV. The 
result agrees with Fig. 1 in Ref. |11| . The peak near zero energy should be 
associated with a small mixture of the spurious mode (pairing rotation). 



In Sees. [2] and El we discuss methods to calculate linear response in nu- 
clei, based on the time-dependent density-functional theory. In this section, 
we will show a feasible real-time method which is, in principle, applicable 
to the non-linear regime as well. 

The time-dependent Hartree-Fock (TDHF) method in the 3D coordi- 
nate representation is a well established method to study nuclear dynamics 
|12| . However, it cannot describe particle-particle (hole-hole) pairing cor- 
relations. The pairing correlations are supposed to be very important not 
only for static properties but also for nuclear dynamics. For instance, it is 
well known that the life time of spontaneous fission is very different between 
even and odd nuclei, which is supposed to be due to the pairing correlations. 
A straightforward extension of the TDHF including the pairing correlations 
is, of course, the time-dependent Hartree-Fock-Bogoliubov (TDHFB) the- 
ory |13| . However it uses the quasi-particle orbitals instead of the occupied 
orbitals whose number is, in principle, infinite. Thus, the accurate calcu- 
lation of TDHFB is presently impractical and a new feasible approach is 
highly desirable. 

In this section, we present the equations of motion of "Canonical-basis 
TDHFB" (Cb-TDHFB) method which we have developed recently [6]. Then, 
we apply the method to the linear-response calculations using the full Skyrme 
functional to show its reliability. For more details, readers should be referred 
to the reference [B]. 



Our starting point is that the TDHFB state can be written in the canon- 
ical form as 



4. Canonical-basis time-dependent HFB method 



4-1. Basic equations 




(12) 



k>0 
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where the creation operator of particles at the canonical state \(f>k(t)) is 
expressed as c\,(t) = J2a I dr<fik(ra;t)ffi (ra). Here, the state k and k are 
not necessarily related to each other by the time reversal, and the time- 
dependent (u, v) factors are complex numbers. Using the density matrix 
and pairing tensor appearing in the HFB equation (pQ), one can write pk(t) = 
\v k (t)\ 2 and K k (t) = u* k (t)v k (t) as 

pk(t) = X^fcMMMOWfcto) = ^{fa{t)\rip^(t){v\fa(t)),(w) 

[IV flV 

K k (t) = Y,^k{t)\p)^- k {t)\y)^ v {t). (14) 

I IV 

Then, utilizing the TDHFB equation, we obtain the following equations for 
the time evolution of pk(t) and n k (t). 

^PfcC*) = Kk(t)AUt)-4(t)A k (t), (15) 
*'^«fe(*) = (Vk(t) + rfc(t)) n k (t) + A k (t) (2p k (t) - 1) , (16) 



where 



A k (t) = -Y,A, u (t)(Mt)\p)(Mt)\"), (17) 

Vk(t) = (Mt)\h(t)\Mt)) + i(^\Mt))- (is) 

So far, there is no approximation in addition to the TDHFB is involved. 
Now, we need to introduce an approximation for the pair potential. Namely, 
the pair potential is assumed to be diagonal in the canonical basis. 

A^{t) = -J2 A k(t)mMt))W\h(t)) - M0*(*)>W*(t)>}- (1 9 ) 

k>0 

In the static limit, this is identical to the BCS approximation. With the 
approximation of Eq. (|19p . one can derive the following simple equations 
for the time-dependent canonical states. 

i^h{t)) = ih{t) - m(t))\Mt)), i^lfaV)) = (Kt) - m(f))\h(t))- 

(20) 

In summary, the Cb-TDHFB equations consists of Eqs. (|20p . (|15p . and 
(|16|) . To derive these equations from the TDHFB equations, we have as- 
sumed the diagonal property of the pair potential, Eq. (fTT)j) , 
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Fig. 4. Calculated (lines) and experimental (symbols) photoabsorption cross section 
in 154 Sm. The SkM* parameter set and the smoothing parameter of V = 2 MeV is 
used. See text for details. 



4-2. Linear response calculation in real time 

We have applied the Cb-TDHFB method to study of electric dipole 
resonances in Ne and Mg isotopes [BJ. Here, we apply the method to GDR in 
the deformed 154 Sm nucleus. We calculate the time evolution of the electric 
dipole moment, starting from the HF+BCS ground state with a perturbative 
instantaneous external dipole field. Then, we perform the Fourier transform 
to obtain the response function. The details of the calculation can be found 
in Ref. [6j. 

We show in Fig. [J]the calculated photoabsorption cross section in 154 Sm. 
Although the pair potential is simplified in the Cb-TDHFB calculation, the 
result is almost identical to the QRPA calculation shown in Fig. [TJ except 
for a small difference seen in the second peak. We have examined the origin 
of this difference and found that the neglect of the residual Coulomb in the 
QRPA calculation is responsible for this small discrepancy. Thus, we may 
conclude that the Cb-TDHFB can reproduce the QRPA result at its small 
amplitude limit. 

It should be emphasized that the computational cost of the Cb-TDHFB 
is significantly smaller than the QRPA. The present calculation in the full 3D 
space can be achieved in roughly 50 CPU hours, while the QRPA calculation 
in Fig. [H that is restricted to the axially symmetric nuclei, requires roughly 
1,000 CPU hours. This is because the Cb-TDHFB treats only the canonical 
states whose number is the same order as the particle number. In contrast, 
in the QRPA (or in the TDHFB), we need to treat the quasiparticle states 
whose number is the same as the dimension of the model space. 
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5. Conclusion 

We have presented our recent developments for studies of nuclear re- 
sponse functions. The parallelized quasiparticle random-phase- approximation 
(QRPA) code is now ready for investigation for heavy axially deformed nu- 
clei. The finite-amplitude method (FAM) was applied to systematic in- 
vestigation of the photoabsorption cross section in light nuclei. Recently, 
the QRPA version of the FAM has been developed for superfluid nuclei, 
including the pairing correlations. We also presented the canonical-basis 
formulation of the TDHFB. This is applicable to large-amplitude nuclear 
dynamics beyond the linear approximation. 
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